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The general relationship between an arbitrary frequency distribution and the expectation value 
of the frequency distributions of its samples is discussed. A wide set of measurable quantities 
( "invariant moments" ) whose expectation value does not in general depend on the size of the sample 
is constructed and illustrated by applying the results to Ewens sampling formula. Invariant moments 
are especially useful in the sampling of systems characterized by the absence of an intrinsic scale. 
Distribution functions that may parametrize the samples of scale-free distributions are considered 
and their invariant expectation values are computed. The conditions under which the scaling limit 
of such distributions may exist are described. 
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. I- INTRODUCTION 



In many interesting physical, biological and social phenomena, whenever no intrinsic scale for the relevant variables 
is present, the emergence of "scaling laws" is phenomenologically observed pj. However, strictly speaking, a power law 
, is not a proper way of fitting empirical data, since no choice of the exponent can keep the higher moments of a power 
law distribution from diverging, while every phenomenological distribution leads to finite values for all moments. This 
is not just a technicality: it is rather a reflection of the fact that the long tail of a power law distribution is in practice 
cutoffed by the existence of some "hidden" scale, irrelevant in the scaling region, but eventually forcing some upper 
limit on the variables describing the system. It would therefore be convenient to parametrize the data by means of 
more regular distribution functions, sufficiently damped for very large values of the variables, but admitting power 
law distributions as regular limits when the control parameter implementing the cutoff is sent to its limiting value. 
A related issue concerns the effects of sampling, which may be non trivial even when we restrict our attention to 
O i' the expectation values of the sampled variables. On average sampling does not affect the distributions of individual 
l— — ', objects belonging to different kinds, but when we consider frequency distributions (that is the number of kinds that 
I ■ are represented k times in a given population) we cannot in general expect that the frequency distribution in the 
samples be the same as in the original population, even after averaging on many different samples, basically because 
£■ — ■ the cutoff induced by sampling acts differently (and in general nontrivially) at different scales. It is therefore quite 
ON ' important to be able to extract from the frequency distribution of the samples some information reflecting directly 
CO some intrinsic property of the underlying distribution. 

Our purpose is therefore threefold. First we want to discuss the general relationship existing between some arbi- 
, trary frequency distribution and the expectation value of the frequency distributions of its samples, and construct 
observables whose expectation values turn out to be independent of the sample size, and therefore coinciding with 
CO | the value taken by the same observables in the full distribution. 

Moreover we want to study classes of distributions whose samples preserve the functional dependence on the 
^ ' parameters present in the original distribution, establishing the connection between the (a priori unknown) values of 
the parameters of the distribution and the (empirically measured) parameters of the sample distributions. 

Finally we want to study the scaling limit of these distributions (when it exists), in order to explore the possibility 
of their use for the phenomenological description of systems that are theoretically expected to show scaling in the 
limit when all empirical cutoffs (including those induced by sampling) are going to disappear. 
In Section |TI] we establish the notation and the general framework of our analysis. 

In Section InTl we construct a wide set of combinations of expectation values that do not depend on the sample size. 
In Section lTVl we apply our approach to the popular Ewens sampling formula, showing that its features are consistent 
with the general pattern and computing its invariant expectation values. 

In Section [V] we consider the limiting case of a small sampling applied to a large population. 

In Section fVTl we focus on the case when the original population and its samples are sufficiently large in comparison 
with typical frequency values, finding a useful mathematical relationship between the generating function of the 
expectation values of the sample distributions and the generating function of the original distribution. 

In Section [VIII we analyze a class of distributions (the so-called negative binomial distributions) admitting a scaling 
limit and enjoying the property that the distribution of expectation values of the samples has the same mathematical 
form as the original distribution. We also compute in a closed form the values ot the basic invariant expectation 
values for these distributions. 
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Finally in Section IVIIII we analyze the scaling limit itself and discuss the conditions under which one may expect 
this limit to be a sensible description of the original system. 

Appendices are devoted to the proofs of some mathematical results and to discussing the issue of correlation between 
random samples. 



II. THE GENERAL FRAMEWORK 



We are considering a set of N objects ("individuals") belonging to S different kinds ("species"), and we assume 
that the set contains N a objects of the a-th kind, subject to the constraint ^ a N a = N ■ 

A sample is a set of n objects, containing h a objects of the a-th kind, subject to the constraint n a = n. 

The probability P{n a y of extracting a specific sample {n a } from a given set {N a } is obtained from the multivariate 
hypergeometric distribution 

We can easily compute the relevant expectation values, obtaining in particular 

(n a ) = N a ^- = np a , (h 2 a ) - (h a ) 2 = — _" np a (l -p a ) 

where p a = N a /N is the probability of extracting an object of the a-th kind in a single extraction. 

It may be useful to consider also the limit of small samples fi a <C N a . In this limit the probability of a specific 
sample is well approximated by the multinomial distribution 

, IL a- 
a—1 

A frequency distribution is a set of values {-/V^}, where is the number of kinds such that for each of them there 
are k objects in the original set. According to the definition, the following conditions must be satisfied: 

N N 



k=l k=l 

The frequency distribution of a sample is a set of values {n{\, satisfying the conditions 

n n 

m = S, ^2 1 n i = n - 

1=0 1 = 1 

Notice that the frequency distribution of a sample formally includes the (unobservable) value no, corresponding to 
the number of kinds, present in the original set, which are not represented in the sample. 

It is in principle possible to compute the probability of any sample distribution {n{\ as a function of a given set 
{Nk}. To this purpose it is convenient to define the intermediate variables Nki, representing the (random) number 
of kinds characterized by k objects in the original set and by I (I < k) objects in the sample. The variables Nki are 
strongly constrained, since they must satisfy all the conditions: 

n N 
1=0 k=l 

The probability P{n m \ of a specific configuration {Nti} follows from the general probability formula 0: 

/Ar\ -1 N r k 1 /A N h - 

subject to the constraint Y^i=o^ kl = 
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The probability P{ ni \ of finding a frequency distribution {ni} in a sample is obtained by summing the probabilities 

P{N k ,\ over all configurations satisfiying the constraint X^feLi ^kl = n i- The corresponding multivariate generating 
function can be defined as 



N k 



£ (™)(M)^^p { „ !} n^ ! - e *wnn* 



{ni} 



1=0 



{N kl } 



k=l 1=0 



It is also possible (and it will be quite convenient) to define a cumulative generating function for the probability of 
finding the frequency distributions P{ ni } for samples of all possible sizes : 
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JV 



^i^peo^w)^ e n(^n^ 



1=0 



xhi 



n[E 

fc=l i=0 



where we used the explicit expression of P{N kt } an d a U the relevant constraints. 

The expectation values (nj) can be computed starting from the above expressions and from the relationship 



N 



N 



k=l k=l{N jm } 



Straightforward manipulations lead to the results [2] 



<JVfci> = N k - 



(»«> 



It is easy to check that the following relationships are satisfied: 

n N n 



N 



k=l 



1=0 



k=l 



In order to fully appreciate the relevance of considerations based on the expectation values we must evaluate the 
weight of the fluctuations. Taking second derivatives of the generating function E(x\ ti) one obtains: 



(nf) - (m) 2 =J2 N kN k 

k,k' 



(N—k—k \ (N—k\ (N—k \ 
\ n-2l J _ \ n-l ) \ n-l ) 

O O (n) 



E^ 



\ n-l ) / fc \ V n-1l I 



Notice that in the large N limit the term quadratic in Nk is depressed by a power of 1/N. This observation suggests 
that very important limits of the above results may be obtained when considering large populations. 



III. INVARIANT EXPECTATION VALUES 



It is very important to be able to define a set of expectation values that are independent of the size of the sample, 
and therefore may reflect very directly the properties of the original frequency distribution. 
Let's consider the following combinations of expectation values: 

where Pi are / arbitrary positive integer numbers, subject only to the constraint that P = '^ i Pi < n. 

The definition of the quantities appearing in the r.h.s. implies that the derivatives with respect to t qi are the joint 
factorial moments of the distribution; therefore we are dealing with weighted combinations of joint factorial moments. 
When some of the indices pi are equal to one, the expectation values may be expressed in terms of a combination of 
lower rank moments (I' (I). 

( p) 

It is possible to recognize that the quantities ( m \ p j) are related (up to a trivial combinatorial factor taking into 
account the existence of n p coincident values of the indices Pj) to the probability of finding the configurations {p{\ in 
the sample containing P elements, and are therefore strictly connected with the probabilities P{n p } ■ 



4 



Exploiting the properties of the generating function E{x\ {ti}) we prove in Appendix [Al that 

for all sets {pi} such that P < n. Hence the expectation values of the nontrivial invariant moments m{ Pi y evaluated for 
samples of arbitrary size n > P, coincide with the moments Ms Pi \ of the original frequency distribution. If the original 
set was generated by a random process, also the M^ pi y will be expectation values. Recalling that e^ N \{ti} = Yl k t^ k 
we may now generate a representation of all P{n p } m terms of Nk, without making use of the coefficients Nki- 

The properties of the binomial coefficients make it possible to invert the relationship between invariant moments 
and joint factorial moments, thus finding that 



I r, I / \ / \ I 

I Ti; \ I 71, \ 



»=1 qi {Pi}i=l V 7 V 7 {p*}*=l 



The basic invariant moments are 

/ \ -i 

E 



m 



According to the inversion formula 



Q £ (n) 



<»<> = Ti-ku-) - Lnrggw - y-'r'Q (>»• 

One may define generating functions for the expectation values of rii and of the basic invariant moments: 

/(«)(*) . VJ (ni )^ ff W (2[) s /(n) (l + £) = ^ (»W£ 

Notice that a special case of the above formula is 

N N /AT\ 

F(t) =Y,N k t\ G( Z ) = F(l + ±) = £ 

fc=0 p=o VF/ 

It is immediate to recognize that g (n) (0) = G(0) = M = S, and ^^(O) = 3£(0) = A/i = 1 



IV. APPLICATION TO EWENS SAMPLING FORMULA 

The multivariate Ewens distribution 0, 0], called in genetics the Ewens sampling formula, describes a specific 
probability for the partition of n into parts, and found its main applications in the context of the neutral theory of 
evolution and in the unified neutral theory of biodiversity [5, 6]. Since Ewens formula and its possibile generalizations 
have been the subject of a wide and still growing literature [7H10| , it may be interesting to apply the results presented 
in Section [IIII to this specific instance. In our notation Ewens probability distribution takes the form 



p M = 



-ft-f-V 



_ r(0 + n) 



where < < oo and ^2 I n i = n - 

Joint factorial moments of the Ewens distribution are easily computed [ll| and one can show that 



n[£]E'wiK-^nG 

1=1 H {n;} (=0 I 



where Q = J2i Qi < 
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We are then left with the task of computing the summations appearing in the equation 



\= ( n \ vnf^ * n ~ Q TT ( 6 \ = p - {n ~ p)lr{9) TT ( t\ V TT ( qi - 1\ r(g + n - Q) 

We prove in Appendix [B] that 



Q)! 



2^ 11 L, _ iy r(6»)(n-Q)! r(6> + p) (n-P)! 



{?i>P;} i=1 

hence 



showing explicitly that the expectation values of the invariant moments of the Ewens distribution are independent of 
the sample size and related to the probability of the configuration {pi} in the sampling of P elements. 

We stress that invariant moments, because of their independence from the size of the sample, may become a highly 
valuable tool in testing the applicability of Ewcns distribution (and of the conceptual assumptions underlying its 
derivation) to the interpretation of actual empirical data. 



V. LARGE POPULATION AND SMALL SAMPLES 



A significant simplification occurs when N — > oo while all other variables are kept finite. Setting x = Nx and to = f 
in the cumulative generating function and taking the large N limit we obtain 



oo ~ n oo ft: 

E(x; M) EE 1 + ]T * n H{U})^ = II I 1 + E 



k\ ( x \ 1 1* 
U 



I J \N 



71=1 k=l 1=1 v 7 k=l 

Let's now define (for n, I different from zero) the following set of coefficients: 

f /f\"i 

71;! 



oo oo 



i=i 



N J 11 



1 1 

({ni} ) S (-i rl ( fl -i)in^y 

i=i 



where s — rii and n = ^ I rii, and notice that the definition of E implies that 

oo n 

\nE(x; {*«}) ee [E c(K) (W}) n( e " ( ° 

n=l { ni } 1 = 1 

It is also possible to recognize that, under the same assumptions, 

oo n 

lnB(£;{t,}) = £[5>< n H{»i})II* 



-- 1 M 



i=i 



where we introduced the large N limit of the basic invariant moments: fhp ee J2 q Nq(q/N) 1 
Comparing the two results we conclude that, for each value of n)0, 

n n 

E c(n) «*}) n( £ ~ (o = [E^^^n^]^- 

{ni} 1 = 1 {rii} 



1=1 



These equations allow in principle for the recursive determination of all e^ n \{ti} in terms of {jrip } (with p < n), 
starting from the initial condition e*- 1 ^ = ti. Higher rank invariant moments (1)1) in the large N limit become 
polynomials in the basic moments. However one must keep in mind that, when the set {N/-} not fixed, but generated 
by a probability distribution (as in the case of the Ewens formula), the expectation values of the products of basic 
moments appearing in the l.h.s. do not coincide with the products of the expectation values. 
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VI. LARGE POPULATIONS AND LARGE SAMPLES 



When k,l <C N,n one may systematically exploit the property that, for small a and b, 

'N - a\ p b {l - p) b ~ a 
n-b) ~^ p n (l - p) N ~ n 

Expressing N and n in terms of N^i one may then obtain 

N r fe 



n 

N' 



p l (l-p) 



k-l 



Pkl = 

k=l L 1=0 

As shown in AppendixfC] when computing expectation values with the above probability distribution, the constraint 
^Ini = n becomes irrelevant in the large N limit, and expectation values of products of Nki with different values 
of the index k factorize into products of expectation values computed for each separate value of k. We can therefore 
compute directly the generating function for a fixed sample size, generalizing the multivariate multinomial distribution: 



N 



1 



\N kl 



k=l{N kl } 



1=0 



N k 

n[5>* 

fc=l 1=0 



A consistency check is easily obtained by observing that e^(l) = 1, because of the property that Yli=o Pkl = 1. 
The expectation value of n\ turns out to be: 



JV 



(ni 



} = j2 N kPki, 



k=l 



and one may check that the conditions on J2i( n i) an d on J2i ' ( n l) are s ^i\\ satisfied. 
We can also estimate the behavior of fluctuations when k.l <S N, n: 



N 



i) - (ni) 2 =J2 N kW ~ Pkl)- 



k=l 

The above expression is always smaller than (m) and as a consequence fluctuations become unimportant for suffi- 
ciently large values of (ni). 

In the same limit we may derive a notable relationship between the generating function of the original frequency 
distribution and the generating function of the expectation values of its samples. In fact we may recognize that for 
sufficiently large N and n 

oo p 

gW(z) = Y,M P ^ = G(z). 

p=0 p ' 

Since in general f^ n \t) = g^(n(t — 1)) and F(t) = G(N(t — 1)), it is then easy to check that in the limit under 
consideration 

&Ht) = F(l-p + pt), 

As a direct consequence of these results, whenever the (size-independent) function 7(2) = G(z) — G(0) can be cast 
into a form exhibiting no explicit parametric dependence on N, the expectation values (ni) can be obtained from Nk 
by the replacement N —¥ n. 

(n) 

Notice that in the limit k,l N,n the definition of the basic invariant moments m p simplifies to 

T 



mW -> P 



l=p 

It is worth analyzing in this limit the explicit expressions of the second basic invariant moment: 

*^Ew-«-E(§)'-i-a^-;.i 

k=l a=l a=l 

As shown in Appendix [D] the above results may be used also in order to parametrize the expected correlation 
between samples under the assumption of independent random sampling. 
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VII. A CLASS OF DISTRIBUTIONS AND ITS PROPERTIES 



As mentioned in the Introduction, distributions found in samples may often correspond to systems whose asymptotic 
(iV — > oo ) distribution is expected to obey a scaling law. However the exponent of the scaling law will in general 
be nontrivial, in contrast with the prediction offered by the simplest neutral models. An example of empirical and 
theoretical evidence for nontriviality is offered by surname frequency distributions (see Ref. [l2j for a recent review), 
recalling that surnames are expected to mimick the behavior of selectively neutral alleles, it is therefore especially 
interesting to consider parametrizations that may reflect notriviality of exponents, ad in particular the class of negative 
binomial distributions [l3j , which can be obtained starting from the generating function 



where 0(a;(l and the parameter c is assumed to vary in the range < c(l. 

The asymptotic behaviour of the distribution for large k is easily obtained by observing that in this limit 

r(fe-c) 1 AT Nil-x) 1 -" x k 



fe! fc 1+c ' x r(i-c) fc 1+c ' 

We can now compute the generating function for the expectation values of the samples according to the general 
rule previously discussed, and obtain 

f c (t) = /c(o) + - (1 ~ y)1 ~ e [i - (1 - yt) c ] , 
y c 

where we have defined y = — —, — . 

a l—x+px 

The distribution of the samples has the same form as the original distribution, once the replacements N — > n and 
x — > y have been performed, and therefore we obtain the asymptotic behaviour 

n(l~y) 1 - c y k 

ni — > 



y r(l-c) /fc 1+c ' 

It is possible to define a combination of parameters independent of the dimension of the sample: 

r at 1 ~ x l ~y 

p = N = n . 

x y 

It is useful to represent x and y in a form showing explicitly their dependence on the dimension of the sample and 
on the invariant parameter /?: 

N n 

y = 



It is now possible to evaluate the expectation value of the invariant moments from the expression 


c 



lc {z) = G c (z) - G c (0) = f [l-(l + |) 



showing no explicit parametric dependence on N; we therefore obtain (for p / 0) 

„_„ r(p- c ) p 1 -* gi l _ /3 



r(i-c) r(i- c )p 1+c ' m 2 i-c' 

The limit of the above results when c — > is smooth, and it corresponds to Fisher distribution fl4| . such that 



and 



F (t) = -/3kL(l-xt), N k = jSy, 



ft 4- N v l 

fo(t) = PHj^)-PHi-yt), ni = Pj- 



The generating function of the invariant moments is obtained from 70(2) = — /31n(l + z//3), and as a consequence 
the expected values of the invariant moments (p ^ 0) are exactly M p = (p — 1)! 

Notice in particular the relationship ft = a, peculiar to Fisher distribution. By comparing these results with the 
large 9 limit of the invariant moments of Ewens distribution we may easily check Watterson's and Hubbell's Q 
observation that in this limit 9 is strictly connected to Fisher's a. 
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VIII. THE SCALING LIMIT 



Let's now consider very large systems, and assume that we can gather information only through the sampling of n 
objects belonging to the system, with n large but not necessarily comparable to N. 

The analysis of the invariant moments may then allow us to check the applicability of a phcnomcnological description 
of the samples based on some distribution falling into the classes discussed in the previous Sections. In the case of a 
positive response to the check it is then possible to find numerical estimates of the parameter j3 and of the exponent 
c. Such estimates are clearly meaningful only if (3 does not turn out to be significantly greater than n. 

Under these assumptions, we can infer a description of the original system, and in the case N))n such a description 
will correspond to computing the limit x — > 1 of the previous results. As a consequence, at least for observable (i.e. 
not too large) values of k, the original distribution is expected to be well described by the scaling form 

N c ^- C 1 

k ^ r(i-c) k 1+c ' 



Appendix A: Expectation values of the invariant moments 



Let's consider the cumulative generating function for the expectation value of a given invariant moment (setting 
P = J2Pi) f° r samples of all possible sizes: 



Let's now observe that, due to the property that 



fe=l i=o ^ ' 



the derivatives appearing in the above defined cumulative generating function, computed at U = 1, can be expressed 
as summations (over {ki} indices) of products of Nk ( times a universal x-dependent factor 



[n 



Q 



(1 + x) 



N-K 



where Q = ^ H an d K — J^i ki- However the following identity holds: 

ki Pi 



Qi 

Pi J \Qt 



PiJ \Qi - Pi 



As a consequence the cumulative generating function is proportional to the factor 



[n 

»=i 



P, 

qi - Pi 



rQ- p 



l-ttl 



aT(l + x) 



N-P 



The summations over the indices {ki} may now be formally performed, and, by matching the coefficient of x N in 
the two sides of the equation, the result can be easily recognized to coincide with the expected value of the invariant 
moment computed for the original distribution times the combinatorial factor (p) . 

In conclusion we find 



N 

E 

n=0 



N\ (n 



n 



(n) \ n 



N 



(m\ N J } )x-(l + x) 



N-P 



Expanding the r.h.s. in powers of x and noticing that ( ) (p) = ( P )( n _ P ) we finally obtain 



N 

E 

n=0 



N\ [N - P 



n - 



p)K;!>^=e( 



N\(N -P 

n- P 



implying immediately that, as long as n > P, 
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Appendix B: Proof of a combinatorial formula 

The identity 

' % - \\ T(6 + n-Q) T(6 + n) 



\j>i-ljT(e)(n-Q)\ T(6 + P)(n-P)V 
where Q — q\ and P — Pi, can be proven by recalling that for real numbers a and positive integers q > p 



1 (a m! "^Vp— 1 



Hence expanding in powers of x the two sides of the identity 



[\{(i~xyv%i-x)- s = {i-xy 



(0 I /') 
i=l 

and exchanging the order of summations in the l.h.s. we get 



r V i ii VP* - V r(0) (m + P - Q)\ ~ T(9 + P)m\ ' 



The desired result is then obtained by setting m = n — P. 

Appendix C: Fluctuations of the sample size in the unconstrained large N limit 

The multivariate generating function e({i;}) was computed in Section [VII after relaxing the constraint Ini = n. 
In order to show that the constraint is automatically satisfied in the large N limit we construct a generating function 
for the expectation value of the powers of v = J2i I n i m the large N distribution: 

im=ev = e p wu}^ ini =n(E ^n^^)^' 

{N kl } fc=l {JV fcl } '=0 

and applying the multinomial formula we obtain 

N k N 

vm = nE^T* = uo--p+p w ) kNk = (o--p+pv) n - 

k=l 1=0 fc=l 

Expanding the result in powers of w we immediately obtain Pi v \ — P/w. 

Since v is distributed according to the binomial distribution, the relevant expectation values are 

i v \ n i y2 \ I v \2 1 \ 

Hence fluctuations of v/N around p vanish like 1/N in the large N limit. 

Appendix D: Correlation between samples 

An important test of randomness in sampling is offered by the measure of the correlation between two different 
samples. Let's consider two random samples, characterized by the sets of values {n a } and {rh a } and by their sizes n 
and m. The index a labels different kinds, as in Section [TTJ The correlation between the two samples is 

c = Ef=i h <^m a 
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Replacing n a and f? a with their expectation values, computed in Section [IIJ we obtain (in the large N limit) 

s s 



a—1 a—1 



a—1 a—1 a—1 a—1 

By making use of the results presented in Section IVII we can now express the expected value of the correlation 
between samples in the form 



(C) 



- + - 

a N 



1 + 1 l + ± 

ct n V a m 



For samples of equal size n the expected value of the correlation takes the form (C) = — 



n a+N 
+n N ' 
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